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ABSTRACT 

The first results from observations of the high mass X-ray binary LS 5039 using the Fermi Gamma-ray 
Space Telescope data between 2008 August and 2009 June are presented. Our results indicate variability that 
is consistent with the binary period, with the emission being modulated with a period of 3.903 ± 0.005 days; 
the first detection of this modulation at GeV energies. The light curve is characterized by a broad peak around 
superior conjunction in agreement with inverse Compton scattering models. The spectrum is represented by 
a power law with an exponential cutoff, yielding an overall flux (100 MeV - 300 GeV) of 4.9 ± 0.5(stat) ± 
1.8(syst) xlO"^ photon cm"^ s"', with a cutoff at 2.1 ± 0.3(stat) ± l.l(syst) GeV and photon index T = 1.9 
± O.l(stat) ± 0.3(syst). The spectrum is observed to vary with orbital phase, specifically between inferior 
and superior conjunction. We suggest that the presence of a cutoff in the spectrum may be indicative of 
magnetospheric emission similar to the emission seen in many pulsars by Fermi. 

Subject headings: binaries: close — gamma rays: observations — stars: variables: other — X-rays: binaries 
— X-rays: individual (LS 5039) 
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1. INTRODUCTION 

LS 5039, PSR B1259-63 and LS I +61°303 are the only 
binaries, with high-mass companions, long known to be spa- 
tially coincident with sources detected at energies greater than 
100 MeV, e.g., those listed in th e Third Energetic G amma-Ray 
Experiment (EGRET) catalog dHartman et al.|[T999.) . The lat- 
ter binary was detected by the Large Area Telescope (LAT) on 
the Fermi mission, confirming it as a GeV gamma-ray source 
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and finding variability that is consi stent with modulation on 
the binary period of 26.6 ± 0.5 days dAbdo et al. 2009b). This 
constituted the first detection of orbital periodicity in high- 
energy (HE) gamma rays (20 MeV - 100 GeV). In this Letter, 
we present the results of Fermi observations of LS 5039. 

LS 5039 is one of a handful of X-ray binaries that have been 
detected recently at very HE 7-rays; results from ^ 70 hr of 
observations distributed over many or bital cycles have been 
presented by lAharonian et al.l (l2005bi [2006). These obser- 
vations yielded a modulation of the very high energy (VHE, 
>100 GeV) gamma -ray flux with a period of 3.9078±0.0015 
days ( Aharonian et al] [20061). consi stent with the orbital pe- 
riod as determined bv lCasares et all 6200 5) from optical spec- 
troscopy. Short timescale variability displayed on top of this 
periodic behavior, both in flux and spectrum, was also re- 
ported. 

The nature of the LS 5039 compact object is unknown: a 
black hole and a neutron star have been invoked as possi- 
ble compact object companions in a slightly eccentric (e ^ 
0.35), 3. 90603 ±0.0 00 17 day orbit around the 06.5V star 
jCasares et al.ll2005h . The discovery of a jet-like radio struc- 
ture in LS 5039 coincident with an X-ray and EGRET source 
prompted a microquasar interpretation (Paredes et al. 2000a). 
Variab ility in the EGRET source could not be established pre- 
cisely ( iTorres et al.ll2001l: iNolan et al.|[2003.) . Recently it has 
been s hown that the X-ray fl ux is modulated on the orbital 
period (Takahashi et al."2009^. 

Ribo et al. (2008) provided Very Long Baseline Array 
(VLBA) radio observations of LS 5039 with morpholog- 
ical and astrometric information at milliarcsecond scales 
that cannot easily be e xplained by a microquasar scenario. 
iMartocchia et al.l (120051) assessed the low X-ray state, show- 
ing the absence of accretion features in the X-ray spectra. 
Thus, measurements at radio and VHE 7-rays in the cases of 
LS I H-6r303 (Dh awan e t al. 2006; Albert et al. 2008) or PSR 
B1259-63 (Ah aronian et a l. 2005a), whose overall spectral 
energy distribution is similar to that of LS 5039, gave the im- 
pression that all three systems are different r ealizations of t he 
same scenario: a pulsar-massive star binary (' Dubus|[2006bb . 

Theoretical computations of the gamma-ray emission in 
both compact object scenarios have been made, with gamma- 
rays produced by inverse Compton (IC) scattering of the 
stellar light by VHE electrons accelerated in the vicinity of 
the compact object. In the case of the black hole com- 
panion, HE and VHE e mission would re s ult from parti- 
cles accelerated in thejet (lBednarekll2007at lBottcheii[2007[: 
Khangulyan et al. 20()8^. Alternatively, it would involve the 
relativistic wind of a young, rotation-powered pulsar, either 
as a result of particle acceleration in the wind interaction re- 
gion ([P ubus 2006b) or b y processes wi thin the pulsar win d 
dSierpowska-Bartosik & ToiTesll2008a[[bl:ICerutti et al.ll2008h . 

2. DATA REDUCTION AND RESULTS 

The LAT onboard Fermi is an electron-positron pair pro- 
duction telescope, featuring solid state silicon trackers and 
cesium iodide calori meters, sensiti ve to photons from ~ 20 
MeV to > 300 GeV (lAtwood et al.l l2009). 

The analysis dataset spanned 2008 August 4, through 2009 
June 22. The data were reduced and analyzed using the Fermi 
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Fig. 1. — Left: the smoothed counts map for 100 MeV-300 GeV of a 10° X 10° region around the LS 5039 location (marked with "x"); the dashed black 
circle indicates the 0.°925 timing analysis aperture. A 0.°3 gaussian smoothing function was applied to the 0.°1 bins. PSR J1826-1256 is marked with "+". 
The exposure varies by less than 7% across the field at a representative energy of 10 GeV. Right: residuals of left panel after subtracting all modelled sources 
excluding LS 5039 and excluding events which arrive during the peaks in the pulsar phase cycle of PSR J1826-1256 (see section|4). The white circles indicate 
the 68% LAT containment region at three energies: 0.3, 1, and 10 GeV. Note that the color scales are different for the two panels. 



Science Tools v9rl5 package^^. The standard onboard filter- 
ing, eve nt reconstruction, an d classification were applied to 
the data ('Atwo od et al.ll20()9h . and for this analysis the high- 
quality ("Pass 6 diffuse") event class is used. Throughout the 
analysis, the "Pass 6 v3 Diffuse" (P6_V3_DIFFUSE) instru- 
ment response functions (IRFs) are applied. 

Time periods when the region around LS 5039 was ob- 
served at a zenith angle greater than 105° and for observa- 
tory rocking angles of greater than 43° were also excluded 
to avoid contamination from the Earth albedo photons. With 
these cuts, a photon count map of a 10° region around the 
binary is shown in Figure [T] 

LS 5039 is detected at a level of 28.5cr (see section |4|i. 
The gtf indsrc tool finds a best-fit location for LS 5039 of 
R.A. = 18^26'"24?7, decl. = -14°48'39"(J2000) with a 95% 
error of 0.054°(including a 20% systematic error derived from 
the internal Fermi 1 1 -month catalog). The nominal position of 
LS 5039 at R.A = 18'' 26'"15.n, Dec. =-14°50'54."3 (J2000) 
(IZacharias et al.l2004ft lies just on the Fermi 95% contour; the 
nominal position was used throughout the analysis. 

3. TIMING ANALYSIS 

LAT light curves were extracted using aperture photome- 
try. The LAT point-spread function (PSF) is strongly energy 
dependent and, particularly since LS 5039 is located in the 
Galactic plane, there is also significant contribution to flux 
within an aperture from diffuse emission and point sources 
that depends on the aperture size and the energy range used. 
The aperture and energy band employed were independently 
chosen to maximize the signal-to-noise level. The optimum 
aperture radius was found to be approximately 0.°925 in the 
energy range 100 MeV-10 GeV. The time binning of the 
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5' See the FSSC website for details 
[http://fermi.gsfc.nasa.gov/ssc/data/analysis/ 



of the Science Tools: 



Fig. 2. — Power spectrum of the light curve. The arrow indicates the known 
orbital period of 3.90603 days i Casares et al. 200^. The dashed lines show 
the 99.99% and 99.99999% significance levels. 

light curve was 1000 s. Exposures were calculated using 
gtexposure assuming a power-law spectrum with a pho- 
ton index of F = 2.5. 

A search was made for periodic modulation by calculat- 
ing the weighted periodogram of the light curve (Lomb l 19761 : 
Scargle 1982; Corbet & Dubois 2007). Since the exposure 
of the time bins was variable, the contribution of each time 
bin to the power spectrum was weighted based on its relative 
exposure. The periodog ram is shown in Figure|2] The arrow 
marks the lCasares et al.l (l2005() orbital period and a highly sig- 
nificant peak is detected at this period; the false alarm prob- 
ability is ^10"'". The significance levels marked are for a 
"blind" search with 5000 independent frequency steps, how- 
ever, the effects of the tuning of the aperture radius and en- 
ergy range are not taken into account. The period error was 
estimated using a Monte Carlo approach; fight curves were 
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simulated using the observed LS 5039 light curve and ran- 
domly shuffling the data points within their errors, assuming 
Gaussian statistics. The corresponding periodogram was then 
calculated and the location of the peak at ^3.9 days recorded. 
From ^-^200,000 simulations, we calculate an error estimate of 
the orbital period of 3.903 ± 0.005 days (Ict). 

4. SPECTRAL ANALYSIS 

The gtlike likelihood fitting tool was used to perform 
the spectral analysis, wherein a spectral-spatial model con- 
taining point and diffuse sources is created and the param- 
eters obtained from a simultaneous maximum likelihood fit 
to the data. The 10° region around the source was modeled 
for Galactic and isotropic diffuse contributions and 22 addi- 
tional significant point sources (taken from the internal Fermi 
11 -month catalog) were included; point sources were mod- 
eled with simple power laws with the exception of the bright, 
nearby pulsar (PSR J 1826- 1256) for which a power law with 
an exponential cut-off was used. The flux contribution of 
PSR J 1826- 1256 to the region was minimized by excluding 
events which arrive during the peaks in the pulsar phase cycle 
( 0.175< <0.3 and 0.625 < <0.775 were excluded; see 
lAbdo et al.ll2009a l). A scaHng factor of 1/0.725 is applied to 
measured fluxes to account for livetime loss due to this cut. 

The 10° region was chosen to capture the broad PSF ob- 
tained at 100 MeV. An alternate fitting method using energy- 
dependent regions of interest was used, yielding compatible 
results that were folded into the systematic errors. 

In analyzing the emission of LS 5039, we used models for 
the Galactic diffuse emission (gll_iem_v02.fit) and isotropic 
backgrounds currently recommended by the LAT team*". The 
model for the Galactic diffuse emission was developed using 
spectral line surveys of H I, CO (as a tracer of H2) to derive 
the distribution of interstellar gas in Galactocentric rings. 

The model of the diffuse gamma-ray emission was then 
constructed by fitting the gamma-ray emissivities of the rings 
in several energy bands to the LAT observations. The fitting 
also require d a model of the IC emission calculated using 
GALPROP (IStrongll2007h and a model for the isoti'opic dif- 
fuse emission. The latter was fitted to the LAT data using an 
analysis of the sky above j/jj = 30 ° and includes the signif- 
icant contribution of residual (misclassified) cosmic rays at 
HEs for the current IRFs. 

Initially a simple power law, , was fitted to the data 
from all orbital phases yielding a photon index of F ~ 2.54. 
However, as clearly seen in Figure [3] the energy spec- 
trum appears to turn over at energies above ~ 2 GeV. The 
possibility of an exponential cutoff was investigated in the 
form £'"'"exp[-(£'/£'cutoff)]- The chance probability to in- 
correctly reject the power-law hypothesis was found to be 
1.6 xlO~'*. The maximum likelihood exponential cutoff 
spect ral model has a likelihood test statistic ( M attox et all 
Il996^ value of ~814.6, or 28.5(t. The photon index is F = 
1.9 ±0.1 (Stat) ± 0.3(syst); flie 100 MeV- 300 GeV flux is 
(4.9 ± 0.5 (stat) ± 1 .8 (syst)) x 10"^ ph cm'^ s"' and the cutoff 
energy is 2.1 ±0.3 (stat) ± 1.1 (syst) GeV (see below for a dis- 
cussion of systematics). The correlation beween the photon 
index and cutoff energy was explored by fitting a family of 
models over a grid of indices and cutoff energies centered on 
the best-fit parameter values. The 1, 2 and 3(t error contours 
are shown in Figured both parameters are bounded and well 

Descriptions of the models are available from the FSSC: 
[http://fermi.gsf c.nasa.gov/ 
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Fig. 3. — Fitted spectrum of LS 5039. Fermi data points are from like- 
lihood fits in each energy bin. The black points (dotted line) represent the 
phase-averaged Fermi/LAI spectrum. The red data points (dotted line) repre- 
sent the spectrum (overall fit) at inferior conjunction (Phase 0.45-0.9); blue 
data points (dotted line) represent the spectrum (overall fit) at superior con- 
junction (Phases, <0.45 and >0.9). Data points above 100 GeV are taken 
from H.E.S.S. observations fAharonian et al. 2006); the data from H.E.S.S. 
are not contemporaneous with Fermi, though they do cover multiple orbital 
periods. 
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Fig. 4. — The Icr, 2cr, and 3a contours for the photon index and cutoff 
energy spectral parameters from fitting to the phase-averaged Fermi data. 

constrained. A total of 359,789 photons were found in the 
10° region. Evaluating the fit parameters, 3992 ± 63 photons 
were observed from LS 5039 above 100 MeV. 

A number of effects are expected to contribute to the sys- 
tematic errors. The largest is uncertainty in the diffuse mod- 
eling as evidenced by the intense swath of photons along the 
Galactic plane shown in Fig [T] A reasonable range of shape 
difference has been explored using the GALPROP model of 
the region. Both models give reasonable residuals maps and 
show differences of 10%, 37%, and 50%, respectively for in- 
dex, cutoff energy and flux. In all diffuse models tested, the 
exponential cutoff model is a significant improvement over 
the power law. 

The impact of systematic uncertainties due to the IRFs are 
estimated by using outlier IRFs that bracket our nominal ones 
in effective area. These are defined by envelopes above and 
below the P6_V3 IRFs by linearly connecting differences of 
(10%, 5%, and 20%) at log(£/MeV) of (2, 2.75, and 4) re- 
spectively. The variation on the index was 15%; the other 
parameter variations were small compared to those due to the 
diffuse modeling. Other effects considered are: fitting tech- 
nique, cuts applied (minimum and maximum energies), but 
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Fig. 5. — Top: flux vs. orbital phase for 0.1-10 GeV. Bottom: vaiiations 
with orbital phase in the hardness ratio of 1-100 to 0.1-1 GeV. 

they are all within the ranges defined by diffuse modeling and 
bracketing effective area variations. 

4.1. Phase resolved spectral analysis 

We also searched for orbital dependence of the spec- 
tral shape, gtlike fits were performed for each phase 
interval of = 0.1 in the same way as for the phase- 
averaged spectral analysis. Figure |5] shows the hardness ra- 
tio, /^i-i()OGev/^().i-iGeV, as a function of orbital phase. Due 
to limited statistics, bin widths of A0 of 0.2 are used for the 
phase interval = 0.3-0.6, while A(/) = 0.1 are employed for 
the other phase intervals. The spectral shape varied such that 
the spectrum is softer around periastron and is harder around 
a pastron. 

lAharonian et al.l (120061) define two broad phase intervals, 
inferior conjunction (0.45 < (j> <0.9) and superior conjunction 
(0.9 < (j) <0.45), both H.E.S.S. spectra being shown in Fig- 
ure |3] Taking the same phase intervals with the LAT data, 
we find a power-law spectrum with F = 2.25 ± 0. 11 at inferior 
conjunction; an energy cutoff was not statistically justified. 
At superior conjunction a power law with an exponential cut- 
off was preferre d with F = 1.91 ± 0.16 and a cutoff energy 
of 1.9±0.5 GeV. lAharonian et al.l 12006) report spectral vari- 
ability in the source between superior and inferior conjunc- 
tion, however, they do not see any indication of long-term 
variability suggesting that it may be reasonable to compare 
non-contemporaneous observations. 

5. DISCUSSION 

The initial assoc iation of LS 5039 wi th the EGRET source 
3EG Jl 824-15 14 (Parede s et aD l2000b) had remained tenta- 
tive due to the large EGRET error circle and the lack of timing 
signatures. The association was bolstered by the discovery of 
point-like, modu lated gamma-ray emission above 100 GeV 
(lAharonian et al . 2005b, 2006). The Fermi observations en- 
abled the detection of an orbital modulation, indicating that 
the binary is also a source of gamma rays above 100 MeV 
The periods determined independently from the ferm/-LAT 
and H.E.S.S. data are self consistent and compatible with the 
binary period o btained from radial v elocity measurements of 
the companion (Aragona et al.l 120091) . This is only t he sec- 
ond high mass X-ray binary after LS I H-61°303 (.Abdo et alj 



l2009bl) to become a confirmed emitter in the HE gamma-ray 
domain. 

The short orbital period means that the compact object 
passes within a stellar ra dius of the surface of the kTi, « 3 eV, 
Ri, sa \QRq companion (ICasares et alJ 12005). Gamma rays 
emitted in the vicinity of the compact object with energies 
above the 30 GeV thr eshold inevitably pair produce with stel- 
lar photons (see e.g. ' Protheroe & Stanev"1987|; iMoskalenkol 
119951: iBottcher & Der mer 2005; Dubus 20061. Emission in 
the Fermi range is largely unaffected by absorption and should 
allow easier identification of the intrinsic spectrum and vari- 
ability of the source. However, it can be affected by cascading 
of higher-energy photons. 

The HE modulation peaks at ~0.0-0.1, close to superior 
conjunction (0 = 0.06) while the VHE modulation peaks close 
to inferior conjunction (0 = 0.72). The phase difference can 
be explained mostly as a result of the competition between 
IC scattering on HE electrons and VHE pair production, as- 
suming the companion star to be close to the compact object. 
The star provides target photons for both processes with the 
radiation density varying by a factor 4 along the eccentric or- 
bit. IC scattering will vary with radiation density but, because 
the source of seed photons is anisotropic, the flux will also 
depen d on the geometry seen by the observer in non- trivial 
ways ( Khangulvan et alj2008l ; ISierpowska-Bartosik & TorresI 
2008b). The emission is enhanced (reduced) when the highly 
relativistic electrons seen by the observer encounter the seed 
photons head-on (rear-on), i.e., at superior (inferior) conjunc- 
tion. Inversely, VHE absorption due to pair production will be 
maximum (minimum) at superior (inferior) conjunction. The 
phases of minimum and maximum flux in Fermi as well as 
the anti-correlation with H.E.S.S. are consistent with these ex- 
pectations, suggesting IC scattering is the dominant radiative 
process above 100 MeV with the additional effect of pair pro- 
duction above 30 GeV fBednarek 2007b; Dubus et al. 20081 
ISierpow ska-Bartosik & Torres 2008a). It is, however, yet un- 
clear whether the IC VHE emission is mainly produced al- 
ready in the pulsar wind zone of the system, given the high 
opacity already found therein for particles accelerated at the 
pulsar site ( Sierpowska-Bartosik & Torres 2008b) or as a re- 
sult of particle acceleration at the shock formed in the wind 
colli sion region (.Dubus 20 06a) or even well beyond the sys- 
tem dBosch-Ramon et all20 08). 

However, the extension of these pictures from the TeV into 
the GeV domain is undermined by the presence of an expo- 
nential cutoff at a few GeV in the Fermi spectrum. A cutoff 
at ^6.3 GeV was also observed in LS I -i-61°303 indicating 
that this may be a common spectral feature in this class of 
source. The companion star in LS I -i-61°303 is a Be star 
with a dense equatorial disk. Passage of the compact object 
through this disk might have explained the exponential cut- 
off: for instance, this would crush a putative pulsar wind neb- 
ula closer to the neutron star, increasing synchrotron losses 
and introducing a strong dep endence with o rbital phase of the 
electron energy distribution (lDubusll2006bl) . But there is no 
such large, systematic contrast in the density of the stellar 
wind from the 06.5V star along the orbit in LS 5039. The 
presence of a similar cutoff in both systems argues against ex- 
planations related to the properties of the orbit or companion 
star. The cutoff seems to require that the radiative process 
or radiating electrons be different between the HE and VHE 
domains, in disagreement with the picture presented above. 

An intringuing possibility is that the emission in the Fermi 
range from both LS 5039 and LS I +6 1 °303 is magnetospheric 
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emission as seen in the dozens of pulsars that have now been 
detected by Fermi. The typical Fermi pulsar emission has a 
hard power-law spectrum with a photon index w 1.5 and an 
exponential cutoff at « 2.5 GeV. In this case, the cutoff energy 
is thought to be set by the balance between acceleration and 
losses to curvature radiation. The emission should in such a 
case be pulsed, although the orbital motion makes it exceed- 
ingly difficult to find in the Fermi data with no prior knowl- 
edge of the spin period. Pulsar wind emission would dom- 
inate in the neighbouring hard X-ray and VHE bands, sup- 
ported by their similar orbital modulations (Hoffmann et al] 
I2OO9I: iTakahashi et alJi2009i) . It is to be noted that there is, 
however, an issue associated with having dominant magneto- 
spheric emission in the HE band: magnetospheric emission 
is usually thought to be due to curvature radiation and has no 
obvious reason to be modulated with the orbital motion (al- 
though magnetospheric emission from pulsars in close bina- 
ries has hardly been modeled). The dense photon environment 
of the binaries may perhaps introduce differences compared 
to gamma-ray magnetospheric emission from isolated pulsars 
(for instance in the pair multiplicity). Hence, the spectrum 
suggests a magnetospheric emission interpretation, which is 
hard to reconcile with the IC scattering interpretation sug- 
gested by the modulation. Future HE and VHE should lead 
to better constraints on the variability of the spectral parame- 



Abdo, A. A., et al. 2009a, Science, 325, 840 

— . 2009b, ApJ, 701, LI 23 

Aharonian, R, et al. 2005a, A&A, 442, 1 

— . 2005b, Science. 309, 746 

— . 2006. A&A, 460, 743 

Albert, J., et al. 2008, ApJ, 684, 1351 

Aragona, C, McSwain, M. V., Grundstrom, E. D., Marsh, A. N., 

Roettenbacher, R. M., Hessler, K. M., Boyajian, T. S., & Ray, R S. 2009, 
ApJ, 698,514 

Atwood, W. B., et al. 2009, ApJ, 697, 1071 

Bednarek, W. 2007a, A&A, 464, 259 

— . 2007b, A&A, 464, 259 

Bosch-Ramon, V., Khangulyan, D., & Aharonian, F. A. 2008, A&A, 489, 
L21 

Bottcher, M. 2007, Astropart. Phys., 27, 278 
Bottcher, M., & Dermer, C. D. 2005, ApJ, 634, L81 

Casares, J., Ribo, M., Ribas, I., Paredes, J. M., Marti, J.. & Herrero, A. 2005. 

MNRAS, 364. 899 
Cerutti, B., Dubus, G., & Henri, G. 2008, A&A, 488, 37 
Corbet, R., & Dubois, R. 2007, in AIP Copnf. Ser. 921, The First GLAST 

Symp., ed. S. Ritz, P Michelson, & C. A. Meegan (Melville, NY:AIP), 

548 

Dhawan, V., Mioduszewski, A., & Rupen, M. 2006, in Proc. VI Microquasar 
Workshop: Microquasars and Beyond, (Como, Italy), ed. T. Belloni, 
PoS(MQW6)052 

Dubus, G. 2006a, A&A, 451, 9 

— . 2006b, A&A, 456, 801 

Dubus, G., Cerutti, B., & Henri, G. 2008, A&A, 477, 691 



ters along the orbit. This will help resolve whether there are 
several components and what their relative amplitudes are. 



The Fermi LAT Collaboration acknowledges generous on- 
going support from a number of agencies and institutes that 
have supported both the development and the operation of 
the LAT as well as scientific data analysis. These include the 
National Aeronautics and Space Administration and the De- 
partment of Energy in the United States, the Commissariat a 
I'Energie Atomique and the Centre National de la Recherche 
Scientifique / Institut National de Physique Nucleaire et de 
Physique des Particules in France, the Agenzia Spaziale Ital- 
iana and the Istituto Nazionale di Fisica Nucleare in Italy, the 
Ministry of Education, Culture, Sports, Science and Technol- 
ogy (MEXT), High Energy Accelerator Research Organiza- 
tion (KEK) and Japan Aerospace Exploration Agency (JAXA) 
in Japan, and the K. A. Wallenberg Foundation, the Swedish 
Research Council and the Swedish National Space Board in 
Sweden. 

Additional support for science analysis during the oper- 
ations phase is gratefully acknowledged from the Spanish 
CSIC and MICINN, the Istituto Nazionale di Astrofisica in 
Italy and the Centre National d'Etudes Spatiales in France. 

Facility: Fermi 



Hartman, R. C, et al. 1999, ApJS, 123, 79 

Hoffmann, A. D., Klochkov, D., Santangelo, A., Horns, D., Segreto, A., 

Staubert, R., & Puehlhofer, G. 2009, A&A, 494, L37 
Khangulyan, D., Aharonian, F, & Bosch-Ramon, V. 2008, MNRAS, 383, 

467 

Lomb. N. R. 1976, Ap&SS, 39, 447 

Martocchia, A., Motch, C, & Negueruela, I. 2005, A&A, 430, 245 
Mattox, J. R., et al. 1996, ApJ, 461, 396 
Moskalenko, I. V. 1995, Space Sci. Rev., 72, 593 

Nolan, P L., Tompkins, W. F, Grenier, I. A., & Michelson, R F. 2003, ApJ, 
597, 615 

Paredes, J. M., Marti, J., Ribo, M., & Massi, M. 2000a, Science, 288, 2340 

— . 2000b, Science, 288, 2340 

Protheroe, R. J., & Stanev, T. 1987, ApJ, 322, 838 

Ribo, M., Paredes, J. M., Moldon, J., Marti, J., & Massi, M. 2008, A&A, 
481, 17 

Scargle, J. D. 1982, ApJ, 263, 835 

Sierpowska-Bartosik, A., & Torres, D. F. 2008a. ApJ, 674, L89 
— . 2008b, Astropart. Phys., 30. 239 
Strong, A. W. 2007, Ap&SS, 309, 35 
Takahashi, T, et al. 2009, ApJ, 697, 592 

ToiTes, D. F., Romero, G. E., Combi, J. A., Benaglia, P., Andernach, H., & 

Punsly, B. 2001, A&A, 370, 468 
Zacharias, N., Urban, S. E., Zacharias, M. I.. Wycoff, G. L., Hall, D. M., 

Monet, D. G., & Rafferty, T. J. 2004, AJ, 127, 3043 



